Halton Sequence Quasi Monte Carlo

numbers=Flatten[Import["~/cernbox/Working/Teaching/Monte Carlo/Lezioni/numbers.txt","Table"]];

Length[numbers]

Histogram[Take[numbers,{1,10000}],{0.01}]

hl=HistogramList[Take[numbers,{1,10000}],{0.01}]
N[StandardDeviation[hl[[2]]]]

StdDevQ[nP_]:=N[StandardDeviation[counts=(HistogramList[Take[numbers,{1,nP}],{0.01}][[2]])]/Mean[counts]]

StdDevRnd[nP_]:=N[StandardDeviation[counts=(HistogramList[RandomReal[{0,1},nP],{0.01}][[2]])]/Mean[counts]]

hl=HistogramList[RandomReal[{0,1},10^4],{0.01}]
N[StandardDeviation[hl[[2]]]]

Table[StdDevQ[nP],{nP,{1000,10000,100000}}]

Table[StdDevRnd[nP],{nP,{1000,10000,100000}}]

{0.6513389472789296`,0.8989331499509894`,1.0249414863810764`}/{3.278334142132975`,9.855521961383532`,28.646927665047514`}

f[x_]:=(x^0.7)*Sqrt[x+1.0/Power[Pi,2.0]]*Sin[1.0/(x+1.0/Power[Pi,2.0])]

Plot[f[x],{x,0,1}]

f[x]//TeXForm

NIntegrate[f[x],{x,0,1}]

MCIntegral[nP_]:=Module[{$numbers},
$numbers=RandomReal[{0,1},nP];
Mean[Table[Total[Map[f,$numbers]/Length[$numbers]],100]]
]

QMCIntegral[nP_]:=Module[{$numbers},
$numbers=Take[numbers,{1,nP}];
Total[Map[f,$numbers]/Length[$numbers]]
]

ra=RandomReal[{0,1},10^6];

Table[A,3]

Map[f,ra]/Length[ra]//Total,Range[4]

QMCres=Table[{nP,QMCIntegral[nP]},{nP,{1000,10000,100000,300000,999000}}]

QMCdiff=Map[{#[[1]],Abs[#[[2]]-0.40442425489145717`]}&,QMCres]

fit=NonlinearModelFit[Map[{Log10[#[[1]]],Log10[#[[2]]]}&,QMCdiff[[1;;-2]]],{a+b*x,b<0},{a,b},x]

Show[ListPlot[fit["Data"]],Plot[fit["BestFit"],{x,3,6}]]

MCres=Table[{nP,MCIntegral[nP]},{nP,{1000,2000,3000,6000,8000,10000}}]

MCdiff=Map[{#[[1]],Abs[#[[2]]-0.40478932457943095`]}&,MCres]

fit=NonlinearModelFit[Map[{Log10[#[[1]]],Log10[#[[2]]]}&,MCdiff[[1;;-2]]],{a+b*x,b<0},{a,b},x]

Show[ListPlot[fit["Data"]],Plot[fit["BestFit"],{x,3,6}]]

NonlinearModelFit[Log10[MCdiff[[1;;-2]]],{a+b*x,b<0},{a,b},x]